#############################################################################################
#  Replication Code for 
#
#  “Value Conflicts Revisited: Muslims, Gender Equality and Gestures of Respect” 
#  British Journal of Political Science
#  Ivarsflaten, E., Helbling, M., Sniderman, P.M., & Traunmüller, R.
#
#  Author: Richard Traunmüller
#  Date: 22.08.23
#  Contact: traunmueller@uni-mannheim.de
#
#############################################################################################

# Load packages
library(psych)
library(foreign)
library(readstata13)
library(survey)
library(stargazer)
library(brant)

# Set working directory
setwd("~/Dropbox/Muslim_Inclusion_Policy_Survey/BJPS_Handshaking/bjps_handshaking_replication_materials/")
      
# Load data
load("bjps_handshaking_replication_data.RData")

# Define function for SE
stderr <- function(x, na.rm=TRUE) {
  if (na.rm) x <- na.omit(x)
  sqrt(var(x)/length(x))
}

###############################################################################################
# Create Variables

dat$pass <- dat$Q63_10==4 & dat$Q73<=2

dat$female <- ifelse(dat$Q65==2, 1, 0)
dat$age<- ifelse(dat$Q67==-99 | dat$Q67>100, NA, dat$Q67)
dat$abi <- ifelse(dat$Q68==-99, NA, dat$Q68)
dat$abi <- ifelse(dat$abi==4, 1, 0)

alpha(cbind(dat$Q17_1, dat$Q17_2, dat$Q17_3, dat$Q17_4))
dat$imcp <-  (dat$Q17_1 + dat$Q17_2 + dat$Q17_3 + dat$Q17_4)/4 

alpha(cbind(dat$Q40, dat$Q41, dat$Q42))
dat$immig <-  (dat$Q40 + dat$Q41 + dat$Q42)/3 

dat$leftright <- dat$Q23

dat$sonntag <- ifelse(dat$Q48==8 | dat$Q48==9, NA, dat$Q48)
round(prop.table(table(dat$sonntag)), 2)

dat$party_l <- ifelse(dat$Q48==4, 1, 0)
dat$party_cl <- ifelse(dat$Q48==2 | dat$Q48==3, 1, 0)
dat$party_cr <- ifelse(dat$Q48==1 | dat$Q48==5, 1, 0)
dat$party_r <- ifelse(dat$Q48==6, 1, 0)

dat$religious <- dat$Q38

alpha(cbind(dat$Q28_1, dat$Q28_2, dat$Q28_3, dat$Q28_4, dat$Q28_5))
dat$christid <-  (dat$Q28_1 + dat$Q28_2 + dat$Q28_3 + dat$Q28_4 + dat$Q28_5)/5 

dat$Q12_hand1.b <- ifelse(dat$Q12_hand1==2, 1, 0)
dat$Q13_hand1.b <- ifelse(dat$Q13_hand1==2, 1, 0)
dat$Q14_hand1.b <- ifelse(dat$Q14_hand1==2, 1, 0)

dat$Q12_politician.b <- ifelse(dat$Q12_politician==2, 1, 0)
dat$Q13_teacher.b <- ifelse(dat$Q13_teacher==2, 1, 0)
dat$Q14_HR.b <- ifelse(dat$Q14_HR==2, 1, 0)

###############################################################################################
# Table A1: Descriptive statistics for respondents who passed both attention checks

dat2 <- dat[dat$pass==T,]  # keep only attentive & serious respondents

cbind(round(apply(cbind(dat2$female, dat2$age, dat2$abi, dat2$imcp, dat2$immig, dat2$leftright, dat2$party_r, dat2$christid), 2, mean, na.rm=T), 2),
  round(apply(cbind(dat2$female, dat2$age, dat2$abi, dat2$imcp, dat2$immig, dat2$leftright, dat2$party_r, dat2$christid), 2, sd, na.rm=T), 2),
  round(apply(cbind(dat2$female, dat2$age, dat2$abi, dat2$imcp, dat2$immig, dat2$leftright, dat2$party_r, dat2$christid), 2, min, na.rm=T), 2),
  round(apply(cbind(dat2$female, dat2$age, dat2$abi, dat2$imcp, dat2$immig, dat2$leftright, dat2$party_r, dat2$christid), 2, max, na.rm=T), 2))


# Table A2: Descriptive statistics for respondents who failed at least one attention check
dat3 <- dat[dat$pass==F,]  # keep only attentive & serious respondents

cbind(round(apply(cbind(dat3$female, dat3$age, dat3$abi, dat3$imcp, dat3$immig, dat3$leftright, dat3$party_r, dat3$christid), 2, mean, na.rm=T), 2),
      round(apply(cbind(dat3$female, dat3$age, dat3$abi, dat3$imcp, dat3$immig, dat3$leftright, dat3$party_r, dat3$christid), 2, sd, na.rm=T), 2),
      round(apply(cbind(dat3$female, dat3$age, dat3$abi, dat3$imcp, dat3$immig, dat3$leftright, dat3$party_r, dat3$christid), 2, min, na.rm=T), 2),
      round(apply(cbind(dat3$female, dat3$age, dat3$abi, dat3$imcp, dat3$immig, dat3$leftright, dat3$party_r, dat3$christid), 2, max, na.rm=T), 2))

# Table A3: Predicting passing both attention checks (OLS)
m.pass <- lm(pass ~ female + age + abi + imcp + immig + leftright + party_r + christid, data=dat)
summary(m.pass)


# Table A4: Balance check of local community center scenario (means)

(table(dat$Q12_hand1.b[is.na(dat$Q12)==F], dat$Q12_politician.b[is.na(dat$Q12)==F]))

rbind(round(aggregate(dat$pass[is.na(dat$Q12)==F], by=list(dat$Q12_hand1.b[is.na(dat$Q12)==F], dat$Q12_politician.b[is.na(dat$Q12)==F]), mean, na.rm=T), 2)$x,
  round(aggregate(dat$female[is.na(dat$Q12)==F], by=list(dat$Q12_hand1.b[is.na(dat$Q12)==F], dat$Q12_politician.b[is.na(dat$Q12)==F]), mean, na.rm=T), 2)$x,
  round(aggregate(dat$age[is.na(dat$Q12)==F], by=list(dat$Q12_hand1.b[is.na(dat$Q12)==F], dat$Q12_politician.b[is.na(dat$Q12)==F]), mean, na.rm=T), 2)$x,
  round(aggregate(dat$abi[is.na(dat$Q12)==F], by=list(dat$Q12_hand1.b[is.na(dat$Q12)==F], dat$Q12_politician.b[is.na(dat$Q12)==F]), mean, na.rm=T), 2)$x,
  round(aggregate(dat$imcp[is.na(dat$Q12)==F], by=list(dat$Q12_hand1.b[is.na(dat$Q12)==F], dat$Q12_politician.b[is.na(dat$Q12)==F]), mean, na.rm=T), 2)$x,
  round(aggregate(dat$immig[is.na(dat$Q12)==F], by=list(dat$Q12_hand1.b[is.na(dat$Q12)==F], dat$Q12_politician.b[is.na(dat$Q12)==F]), mean, na.rm=T), 2)$x,
  round(aggregate(dat$leftright[is.na(dat$Q12)==F], by=list(dat$Q12_hand1.b[is.na(dat$Q12)==F], dat$Q12_politician.b[is.na(dat$Q12)==F]), mean, na.rm=T), 2)$x,
  round(aggregate(dat$party_r[is.na(dat$Q12)==F], by=list(dat$Q12_hand1.b[is.na(dat$Q12)==F], dat$Q12_politician.b[is.na(dat$Q12)==F]), mean, na.rm=T), 2)$x,
  round(aggregate(dat$christid[is.na(dat$Q12)==F], by=list(dat$Q12_hand1.b[is.na(dat$Q12)==F], dat$Q12_politician.b[is.na(dat$Q12)==F]), mean, na.rm=T), 2)$x)

# Table A5: Balance check of local school scenario (means)

(table(dat$Q13_hand1.b[is.na(dat$Q13)==F], dat$Q13_teacher.b[is.na(dat$Q13)==F]))

rbind(round(aggregate(dat$pass[is.na(dat$Q13)==F], by=list(dat$Q13_hand1.b[is.na(dat$Q13)==F], dat$Q13_teacher.b[is.na(dat$Q13)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$female[is.na(dat$Q13)==F], by=list(dat$Q13_hand1.b[is.na(dat$Q13)==F], dat$Q13_teacher.b[is.na(dat$Q13)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$age[is.na(dat$Q13)==F], by=list(dat$Q13_hand1.b[is.na(dat$Q13)==F], dat$Q13_teacher.b[is.na(dat$Q13)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$abi[is.na(dat$Q13)==F], by=list(dat$Q13_hand1.b[is.na(dat$Q13)==F], dat$Q13_teacher.b[is.na(dat$Q13)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$imcp[is.na(dat$Q13)==F], by=list(dat$Q13_hand1.b[is.na(dat$Q13)==F], dat$Q13_teacher.b[is.na(dat$Q13)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$immig[is.na(dat$Q13)==F], by=list(dat$Q13_hand1.b[is.na(dat$Q13)==F], dat$Q13_teacher.b[is.na(dat$Q13)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$leftright[is.na(dat$Q13)==F], by=list(dat$Q13_hand1.b[is.na(dat$Q13)==F], dat$Q13_teacher.b[is.na(dat$Q13)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$party_r[is.na(dat$Q13)==F], by=list(dat$Q13_hand1.b[is.na(dat$Q13)==F], dat$Q13_teacher.b[is.na(dat$Q13)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$christid[is.na(dat$Q13)==F], by=list(dat$Q13_hand1.b[is.na(dat$Q13)==F], dat$Q13_teacher.b[is.na(dat$Q13)==F]), mean, na.rm=T), 2)$x)


# A6: Balance check of local company scenario (means)

(table(dat$Q14_hand1.b[is.na(dat$Q14)==F], dat$Q14_HR.b[is.na(dat$Q14)==F]))

rbind(round(aggregate(dat$pass[is.na(dat$Q14)==F], by=list(dat$Q14_hand1.b[is.na(dat$Q14)==F], dat$Q14_HR.b[is.na(dat$Q14)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$female[is.na(dat$Q14)==F], by=list(dat$Q14_hand1.b[is.na(dat$Q14)==F], dat$Q14_HR.b[is.na(dat$Q14)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$age[is.na(dat$Q14)==F], by=list(dat$Q14_hand1.b[is.na(dat$Q14)==F], dat$Q14_HR.b[is.na(dat$Q14)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$abi[is.na(dat$Q14)==F], by=list(dat$Q14_hand1.b[is.na(dat$Q14)==F], dat$Q14_HR.b[is.na(dat$Q14)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$imcp[is.na(dat$Q14)==F], by=list(dat$Q14_hand1.b[is.na(dat$Q14)==F], dat$Q14_HR.b[is.na(dat$Q14)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$immig[is.na(dat$Q14)==F], by=list(dat$Q14_hand1.b[is.na(dat$Q14)==F], dat$Q14_HR.b[is.na(dat$Q14)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$leftright[is.na(dat$Q14)==F], by=list(dat$Q14_hand1.b[is.na(dat$Q14)==F], dat$Q14_HR.b[is.na(dat$Q14)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$party_r[is.na(dat$Q14)==F], by=list(dat$Q14_hand1.b[is.na(dat$Q14)==F], dat$Q14_HR.b[is.na(dat$Q14)==F]), mean, na.rm=T), 2)$x,
      round(aggregate(dat$christid[is.na(dat$Q14)==F], by=list(dat$Q14_hand1.b[is.na(dat$Q14)==F], dat$Q14_HR.b[is.na(dat$Q14)==F]), mean, na.rm=T), 2)$x)

###############################################################################################
# Table A7: OLS Results: Should be Asked to Shake Hands [%]. Full Sample

m.1 <- lm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b, data=dat)
m.2 <- lm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b, data=dat)
m.3 <- lm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b, data=dat)

m.1.b <- lm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b + female + age + abi, data=dat)
m.2.b <- lm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b + female + age + abi, data=dat)
m.3.b <- lm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b + female + age + abi, data=dat)

stargazer(m.1, m.1.b, m.2, m.2.b, m.3, m.3.b, type="text", omit.stat=c("adj.rsq", "f", "ser"))

# Table A8: OLS Results: Should be Asked to Shake Hands [%]. Passed Attention Check.

m.1 <- lm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b, data=dat2)
m.2 <- lm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b, data=dat2)
m.3 <- lm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b, data=dat2)

m.1.b <- lm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b + female + age + abi, data=dat2)
m.2.b <- lm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b + female + age + abi, data=dat2)
m.3.b <- lm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b + female + age + abi, data=dat2)

stargazer(m.1, m.1.b, m.2, m.2.b, m.3, m.3.b, type="text", omit.stat=c("adj.rsq", "f", "ser"))

#############################################################################################
# Figure 2: Results of Handshaking Experiments. 

png("bjps_handshaking_figure_2.png", width=3, height=9, units="in", res=600)
par(mfrow=c(3,1), las=1, tck=-.02, mar=c(2, 3, 3, 2))

# Politicians 
exp.1 <- aggregate(dat$Q12>=3, list(dat$Q12_politician, dat$Q12_hand1), mean, na.rm=T)
se.exp.1 <- aggregate(dat$Q12>=3, list(dat$Q12_politician, dat$Q12_hand1), stderr, na.rm=T)

par(mgp=c(2, .5, 0), family="Gill Sans", font.main=1, cex.main=1, cex.lab=1, cex.axis=1)
plot(exp.1[c(1,3),2], exp.1[c(1,3),3], pch=19, type="o", ylim=c(0, 1), xlim=c(0.7, 2.3), axes=F, ylab="Should Be Asked to Shake Hands [%]", xlab="", main="Local Community Center")
points(exp.1[c(2,4),2], exp.1[c(2,4),3], type="o")
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(1, at=c(1, 2), labels=c("No Mention", "Hand on Heart"))
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
segments(1, exp.1[1,3]-1.96*se.exp.1[1,3], 1, exp.1[1,3]+1.96*se.exp.1[1,3])
segments(1, exp.1[2,3]-1.96*se.exp.1[2,3], 1, exp.1[2,3]+1.96*se.exp.1[2,3])
segments(2, exp.1[3,3]-1.96*se.exp.1[3,3], 2, exp.1[3,3]+1.96*se.exp.1[3,3])
segments(2, exp.1[4,3]-1.96*se.exp.1[4,3], 2, exp.1[4,3]+1.96*se.exp.1[4,3])
points(exp.1[c(2,4),2], exp.1[c(2,4),3], pch=19, col="white", cex=.7)
text(1, exp.1[1,3]-1.96*se.exp.1[1,3], "Politicians", pos=1, cex=.8)
text(1, exp.1[2,3]+1.96*se.exp.1[2,3], "Female Politicians", pos=3, cex=.8)

text(1, exp.1[1,3], round(exp.1[1,3]*100), pos=2, cex=.8)
text(1, exp.1[2,3], round(exp.1[2,3]*100), pos=2, cex=.8)
text(2, exp.1[3,3], paste(round(exp.1[3,3]*100), "+"), pos=4, cex=.8)
text(2, exp.1[4,3], paste(round(exp.1[4,3]*100), "**"), pos=4, cex=.8)

box()

# Teachers
exp.2 <- aggregate(dat$Q13>=3, list(dat$Q13_teacher, dat$Q13_hand1), mean, na.rm=T)
se.exp.2 <- aggregate(dat$Q13>=3, list(dat$Q13_teacher, dat$Q13_hand1), stderr, na.rm=T)

par(mgp=c(2, .5, 0), family="Gill Sans", font.main=1, cex.main=1)
plot(exp.2[c(1,3),2], exp.2[c(1,3),3], pch=19, type="o", ylim=c(0, 1), xlim=c(0.7, 2.3), axes=F, ylab="Should Be Asked to Shake Hands [%]", xlab="", main="Local School")
points(exp.2[c(2,4),2], exp.2[c(2,4),3], type="o")
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(1, at=c(1, 2), labels=c("No Mention", "Hand on Heart"))
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
segments(1, exp.2[1,3]-1.96*se.exp.2[1,3], 1, exp.2[1,3]+1.96*se.exp.2[1,3])
segments(1, exp.2[2,3]-1.96*se.exp.2[2,3], 1, exp.2[2,3]+1.96*se.exp.2[2,3])
segments(2, exp.2[3,3]-1.96*se.exp.2[3,3], 2, exp.2[3,3]+1.96*se.exp.2[3,3])
segments(2, exp.2[4,3]-1.96*se.exp.2[4,3], 2, exp.2[4,3]+1.96*se.exp.2[4,3])
points(exp.2[c(2,4),2], exp.2[c(2,4),3], pch=19, col="white", cex=.7)
text(1, exp.2[1,3]-1.96*se.exp.2[1,3], "Teachers", pos=1, cex=.8)
text(1, exp.2[2,3]+1.96*se.exp.2[2,3], "Female Teachers", pos=3, cex=.8)

text(1, exp.2[1,3], round(exp.2[1,3]*100), pos=2, cex=.8)
text(1, exp.2[2,3], round(exp.2[2,3]*100), pos=2, cex=.8)
text(2, exp.2[3,3], paste(round(exp.2[3,3]*100), "*"), pos=4, cex=.8)
text(2, exp.2[4,3], paste(round(exp.2[4,3]*100), "**"), pos=4, cex=.8)

box()


# HR Managers
exp.3 <- aggregate(dat$Q14>=3, list(dat$Q14_HR, dat$Q14_hand1), mean, na.rm=T)
se.exp.3 <- aggregate(dat$Q14>=3, list(dat$Q14_HR, dat$Q14_hand1), stderr, na.rm=T)

par(mgp=c(2, .5, 0), family="Gill Sans", font.main=1, cex.main=1)
plot(exp.3[c(1,3),2], exp.3[c(1,3),3], pch=19, type="o", ylim=c(0, 1), xlim=c(0.7, 2.3), axes=F, ylab="Should Be Asked to Shake Hands [%]", xlab="", main="Local Company")
points(exp.3[c(2,4),2], exp.3[c(2,4),3], type="o")
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(1, at=c(1, 2), labels=c("No Mention", "Hand on Heart"))
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
segments(1, exp.3[1,3]-1.96*se.exp.3[1,3], 1, exp.3[1,3]+1.96*se.exp.3[1,3])
segments(1, exp.3[2,3]-1.96*se.exp.3[2,3], 1, exp.3[2,3]+1.96*se.exp.3[2,3])
segments(2, exp.3[3,3]-1.96*se.exp.3[3,3], 2, exp.3[3,3]+1.96*se.exp.3[3,3])
segments(2, exp.3[4,3]-1.96*se.exp.3[4,3], 2, exp.3[4,3]+1.96*se.exp.3[4,3])
points(exp.3[c(2,4),2], exp.3[c(2,4),3], pch=19, col="white", cex=.7)
text(1, exp.3[1,3]-1.96*se.exp.3[1,3], "HR Managers", pos=1, cex=.8)
text(1, exp.3[2,3]+1.96*se.exp.3[2,3], "Female HR Managers", pos=3, cex=.8)

text(1, exp.3[1,3], round(exp.3[1,3]*100), pos=2, cex=.8)
text(1, exp.3[2,3], round(exp.3[2,3]*100), pos=2, cex=.8)
text(2, exp.3[3,3], paste(round(exp.3[3,3]*100), "**"), pos=4, cex=.8)
text(2, exp.3[4,3], paste(round(exp.3[4,3]*100), "**"), pos=4, cex=.8)

box()
dev.off()

t.test(Q12>=3 ~ Q12_hand1, data=dat[dat$Q12_politician==1,])
t.test(Q12>=3 ~ Q12_hand1, data=dat[dat$Q12_politician==2,])

t.test(Q13>=3 ~ Q13_hand1, data=dat[dat$Q13_teacher==1,])
t.test(Q13>=3 ~ Q13_hand1, data=dat[dat$Q13_teacher==2,])

t.test(Q14>=3 ~ Q14_hand1, data=dat[dat$Q14_HR==1,])
t.test(Q14>=3 ~ Q14_hand1, data=dat[dat$Q14_HR==2,])

##############################################################################################
# Figure A1: Results of Handshaking Experiments by Respondents’ Gender. 

png("bjps_handshaking_figure_A1.png", width=5, height=9, units="in", res=600)
par(mfrow=c(3,1))

# Politicians 
exp.1.f <- aggregate(dat$Q12[dat$female==1]>=3, list(dat$Q12_politician[dat$female==1], dat$Q12_hand1[dat$female==1]), mean, na.rm=T)
se.exp.1.f <- aggregate(dat$Q12[dat$female==1]>=3, list(dat$Q12_politician[dat$female==1], dat$Q12_hand1[dat$female==1]), stderr, na.rm=T)
exp.1.m <- aggregate(dat$Q12[dat$female==0]>=3, list(dat$Q12_politician[dat$female==0], dat$Q12_hand1[dat$female==0]), mean, na.rm=T)
se.exp.1.m <- aggregate(dat$Q12[dat$female==0]>=3, list(dat$Q12_politician[dat$female==0], dat$Q12_hand1[dat$female==0]), stderr, na.rm=T)

par(mgp=c(2, .5, 0), family="Gill Sans", font.main=1, cex.main=1, cex.lab=1, cex.axis=1)
plot(exp.1.f[c(1,3),2], exp.1.f[c(1,3),3], pch=19, type="o", ylim=c(0, 1), xlim=c(0.7, 4.3), axes=F, ylab="Should Be Asked to Shake Hands [%]", xlab="", main="Local Community Center")
points(exp.1.f[c(2,4),2], exp.1.f[c(2,4),3], type="o")
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(1, at=c(1:4), labels=c("No Mention", "Hand on Heart", "No Mention", "Hand on Heart"))
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
segments(1, exp.1.f[1,3]-1.96*se.exp.1.f[1,3], 1, exp.1.f[1,3]+1.96*se.exp.1.f[1,3])
segments(1, exp.1.f[2,3]-1.96*se.exp.1.f[2,3], 1, exp.1.f[2,3]+1.96*se.exp.1.f[2,3])
segments(2, exp.1.f[3,3]-1.96*se.exp.1.f[3,3], 2, exp.1.f[3,3]+1.96*se.exp.1.f[3,3])
segments(2, exp.1.f[4,3]-1.96*se.exp.1.f[4,3], 2, exp.1.f[4,3]+1.96*se.exp.1.f[4,3])
points(exp.1.f[c(2,4),2], exp.1.f[c(2,4),3], pch=19, col="white", cex=.7)
text(1, exp.1.f[1,3]-1.96*se.exp.1.f[1,3], "Politicians", pos=1, cex=.8)
text(1, exp.1.f[2,3]+1.96*se.exp.1.f[2,3], "Female Politicians", pos=3, cex=.8)
text(1, exp.1.f[1,3], round(exp.1.f[1,3]*100), pos=2, cex=.8)
text(1, exp.1.f[2,3], round(exp.1.f[2,3]*100), pos=2, cex=.8)
text(2, exp.1.f[3,3], paste(round(exp.1.f[3,3]*100), "*"), pos=4, cex=.8)
text(2, exp.1.f[4,3], paste(round(exp.1.f[4,3]*100), "**"), pos=4, cex=.8)

points(exp.1.m[c(1,3),2]+2, exp.1.m[c(1,3),3], type="o", pch=19)
points(exp.1.m[c(2,4),2]+2, exp.1.m[c(2,4),3], type="o")
segments(1+2, exp.1.m[1,3]-1.96*se.exp.1.m[1,3], 1+2, exp.1.m[1,3]+1.96*se.exp.1.m[1,3])
segments(1+2, exp.1.m[2,3]-1.96*se.exp.1.m[2,3], 1+2, exp.1.m[2,3]+1.96*se.exp.1.m[2,3])
segments(2+2, exp.1.m[3,3]-1.96*se.exp.1.m[3,3], 2+2, exp.1.m[3,3]+1.96*se.exp.1.m[3,3])
segments(2+2, exp.1.m[4,3]-1.96*se.exp.1.m[4,3], 2+2, exp.1.m[4,3]+1.96*se.exp.1.m[4,3])
points(exp.1.m[c(2,4),2]+2, exp.1.m[c(2,4),3], pch=19, col="white", cex=.7)
text(1+2, exp.1.m[1,3], round(exp.1.m[1,3]*100), pos=2, cex=.8)
text(1+2, exp.1.m[2,3], round(exp.1.m[2,3]*100), pos=2, cex=.8)
text(2+2, exp.1.m[3,3], paste(round(exp.1.m[3,3]*100), "n.s."), pos=4, cex=.8)
text(2+2, exp.1.m[4,3], paste(round(exp.1.m[4,3]*100), "**"), pos=4, cex=.8)

text(1.5, 0, "Female", pos=3)
text(3.5, 0, "Male", pos=3)

box()

t.test(Q12>=3 ~ Q12_hand1, data=dat[dat$Q12_politician==1 & dat$female==1,])
t.test(Q12>=3 ~ Q12_hand1, data=dat[dat$Q12_politician==2 & dat$female==1,])
t.test(Q12>=3 ~ Q12_hand1, data=dat[dat$Q12_politician==1 & dat$female==0,])
t.test(Q12>=3 ~ Q12_hand1, data=dat[dat$Q12_politician==2 & dat$female==0,])

summary(lm(Q12>=3 ~ Q12_hand1.b*female, data=dat[dat$Q12_politician==1,]))
summary(lm(Q12>=3 ~ Q12_hand1.b*female, data=dat[dat$Q12_politician==2,]))

summary(lm(Q12>=3 ~ Q12_hand1.b*dat$Q12_politician.b*female, data=dat))

# Teachers
exp.2.f <- aggregate(dat$Q13[dat$female==1]>=3, list(dat$Q13_teacher[dat$female==1], dat$Q13_hand1[dat$female==1]), mean, na.rm=T)
se.exp.2.f <- aggregate(dat$Q13[dat$female==1]>=3, list(dat$Q13_teacher[dat$female==1], dat$Q13_hand1[dat$female==1]), stderr, na.rm=T)
exp.2.m <- aggregate(dat$Q13[dat$female==0]>=3, list(dat$Q13_teacher[dat$female==0], dat$Q13_hand1[dat$female==0]), mean, na.rm=T)
se.exp.2.m <- aggregate(dat$Q13[dat$female==0]>=3, list(dat$Q13_teacher[dat$female==0], dat$Q13_hand1[dat$female==0]), stderr, na.rm=T)

par(mgp=c(2, .5, 0), family="Gill Sans", font.main=1, cex.main=1)
plot(exp.2.f[c(1,3),2], exp.2.f[c(1,3),3], pch=19, type="o", ylim=c(0, 1), xlim=c(0.7, 4.3), axes=F, ylab="Should Be Asked to Shake Hands [%]", xlab="", main="Local School")
points(exp.2.f[c(2,4),2], exp.2.f[c(2,4),3], type="o")
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(1, at=c(1, 2, 3, 4), labels=c("No Mention", "Hand on Heart", "No Mention", "Hand on Heart"))
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
segments(1, exp.2.f[1,3]-1.96*se.exp.2.f[1,3], 1, exp.2.f[1,3]+1.96*se.exp.2.f[1,3])
segments(1, exp.2.f[2,3]-1.96*se.exp.2.f[2,3], 1, exp.2.f[2,3]+1.96*se.exp.2.f[2,3])
segments(2, exp.2.f[3,3]-1.96*se.exp.2.f[3,3], 2, exp.2.f[3,3]+1.96*se.exp.2.f[3,3])
segments(2, exp.2.f[4,3]-1.96*se.exp.2.f[4,3], 2, exp.2.f[4,3]+1.96*se.exp.2.f[4,3])
points(exp.2.f[c(2,4),2], exp.2.f[c(2,4),3], pch=19, col="white", cex=.7)
text(1, exp.2.f[1,3]-1.96*se.exp.2.f[1,3], "Teachers", pos=1, cex=.8)
text(1, exp.2.f[2,3]+1.96*se.exp.2.f[2,3], "Female Teachers", pos=3, cex=.8)
text(1, exp.2.f[1,3], round(exp.2.f[1,3]*100), pos=2, cex=.8)
text(1, exp.2.f[2,3], round(exp.2.f[2,3]*100), pos=2, cex=.8)
text(2, exp.2.f[3,3], paste(round(exp.2.f[3,3]*100), "n.s."), pos=4, cex=.8)
text(2, exp.2.f[4,3], paste(round(exp.2.f[4,3]*100), "*"), pos=4, cex=.8)

points(exp.2.m[c(1,3),2]+2, exp.2.m[c(1,3),3], type="o", pch=19)
points(exp.2.m[c(2,4),2]+2, exp.2.m[c(2,4),3], type="o")
segments(1+2, exp.2.m[1,3]-1.96*se.exp.2.m[1,3], 1+2, exp.2.m[1,3]+1.96*se.exp.2.m[1,3])
segments(1+2, exp.2.m[2,3]-1.96*se.exp.2.m[2,3], 1+2, exp.2.m[2,3]+1.96*se.exp.2.m[2,3])
segments(2+2, exp.2.m[3,3]-1.96*se.exp.2.m[3,3], 2+2, exp.2.m[3,3]+1.96*se.exp.2.m[3,3])
segments(2+2, exp.2.m[4,3]-1.96*se.exp.2.m[4,3], 2+2, exp.2.m[4,3]+1.96*se.exp.2.m[4,3])
points(exp.2.m[c(2,4),2]+2, exp.2.m[c(2,4),3], pch=19, col="white", cex=.7)
text(1+2, exp.2.m[1,3], round(exp.2.m[1,3]*100), pos=2, cex=.8)
text(1+2, exp.2.m[2,3], round(exp.2.m[2,3]*100), pos=2, cex=.8)
text(2+2, exp.2.m[3,3], paste(round(exp.2.m[3,3]*100), "*"), pos=4, cex=.8)
text(2+2, exp.2.m[4,3], paste(round(exp.2.m[4,3]*100), "**"), pos=4, cex=.8)

text(1.5, 0, "Female", pos=3)
text(3.5, 0, "Male", pos=3)

box()


t.test(Q13>=3 ~ Q13_hand1, data=dat[dat$Q13_teacher==1 & dat$female==1,])
t.test(Q13>=3 ~ Q13_hand1, data=dat[dat$Q13_teacher==2 & dat$female==1,])
t.test(Q13>=3 ~ Q13_hand1, data=dat[dat$Q13_teacher==1 & dat$female==0,])
t.test(Q13>=3 ~ Q13_hand1, data=dat[dat$Q13_teacher==2 & dat$female==0,])

dat$Q13_hand1.b <- ifelse(dat$Q13_hand1==2, 1, 0)
summary(lm(Q13>=3 ~ Q13_hand1.b*female, data=dat[dat$Q13_teacher==1,]))
summary(lm(Q13>=3 ~ Q13_hand1.b*female, data=dat[dat$Q13_teacher==2,]))

summary(lm(Q13>=3 ~ Q13_hand1.b*dat$Q13_teacher.b*female, data=dat))

# HR Managers
exp.3.f <- aggregate(dat$Q14[dat$female==1]>=3, list(dat$Q14_HR[dat$female==1], dat$Q14_hand1[dat$female==1]), mean, na.rm=T)
se.exp.3.f <- aggregate(dat$Q14[dat$female==1]>=3, list(dat$Q14_HR[dat$female==1], dat$Q14_hand1[dat$female==1]), stderr, na.rm=T)
exp.3.m <- aggregate(dat$Q14[dat$female==0]>=3, list(dat$Q14_HR[dat$female==0], dat$Q14_hand1[dat$female==0]), mean, na.rm=T)
se.exp.3.m <- aggregate(dat$Q14[dat$female==0]>=3, list(dat$Q14_HR[dat$female==0], dat$Q14_hand1[dat$female==0]), stderr, na.rm=T)

par(mgp=c(2, .5, 0), family="Gill Sans", font.main=1, cex.main=1)
plot(exp.3.f[c(1,3),2], exp.3.f[c(1,3),3], pch=19, type="o", ylim=c(0, 1), xlim=c(0.7, 4.3), axes=F, ylab="Should Be Asked to Shake Hands [%]", xlab="", main="Local Company")
points(exp.3.f[c(2,4),2], exp.3.f[c(2,4),3], type="o")
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(1, at=c(1, 2, 3, 4), labels=c("No Mention", "Hand on Heart", "No Mention", "Hand on Heart"))
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
segments(1, exp.3.f[1,3]-1.96*se.exp.3.f[1,3], 1, exp.3.f[1,3]+1.96*se.exp.3.f[1,3])
segments(1, exp.3.f[2,3]-1.96*se.exp.3.f[2,3], 1, exp.3.f[2,3]+1.96*se.exp.3.f[2,3])
segments(2, exp.3.f[3,3]-1.96*se.exp.3.f[3,3], 2, exp.3.f[3,3]+1.96*se.exp.3.f[3,3])
segments(2, exp.3.f[4,3]-1.96*se.exp.3.f[4,3], 2, exp.3.f[4,3]+1.96*se.exp.3.f[4,3])
points(exp.3.f[c(2,4),2], exp.3.f[c(2,4),3], pch=19, col="white", cex=.7)
text(1, exp.3.f[1,3]-1.96*se.exp.3.f[1,3], "HR Managers", pos=1, cex=.8)
text(1, exp.3.f[2,3]+1.96*se.exp.3.f[2,3], "Female HR Managers", pos=3, cex=.8)
text(1, exp.3.f[1,3], round(exp.3.f[1,3]*100), pos=2, cex=.8)
text(1, exp.3.f[2,3], round(exp.3.f[2,3]*100), pos=2, cex=.8)
text(2, exp.3.f[3,3], paste(round(exp.3.f[3,3]*100), "**"), pos=4, cex=.8)
text(2, exp.3.f[4,3], paste(round(exp.3.f[4,3]*100), "**"), pos=4, cex=.8)

points(exp.3.m[c(1,3),2]+2, exp.3.m[c(1,3),3], type="o", pch=19)
points(exp.3.m[c(2,4),2]+2, exp.3.m[c(2,4),3], type="o")
segments(1+2, exp.3.m[1,3]-1.96*se.exp.3.m[1,3], 1+2, exp.3.m[1,3]+1.96*se.exp.3.m[1,3])
segments(1+2, exp.3.m[2,3]-1.96*se.exp.3.m[2,3], 1+2, exp.3.m[2,3]+1.96*se.exp.3.m[2,3])
segments(2+2, exp.3.m[3,3]-1.96*se.exp.3.m[3,3], 2+2, exp.3.m[3,3]+1.96*se.exp.3.m[3,3])
segments(2+2, exp.3.m[4,3]-1.96*se.exp.3.m[4,3], 2+2, exp.3.m[4,3]+1.96*se.exp.3.m[4,3])
points(exp.3.m[c(2,4),2]+2, exp.3.m[c(2,4),3], pch=19, col="white", cex=.7)
text(1+2, exp.3.m[1,3], round(exp.3.m[1,3]*100), pos=2, cex=.8)
text(1+2, exp.3.m[2,3], round(exp.3.m[2,3]*100), pos=2, cex=.8)
text(2+2, exp.3.m[3,3], paste(round(exp.3.m[3,3]*100), "*"), pos=4, cex=.8)
text(2+2, exp.3.m[4,3], paste(round(exp.3.m[4,3]*100), "*"), pos=4, cex=.8)

text(1.5, 0, "Female", pos=3)
text(3.5, 0, "Male", pos=3)

box()

t.test(Q14>=3 ~ Q14_hand1, data=dat[dat$Q14_HR==1 & dat$female==1,])
t.test(Q14>=3 ~ Q14_hand1, data=dat[dat$Q14_HR==2 & dat$female==1,])
t.test(Q14>=3 ~ Q14_hand1, data=dat[dat$Q14_HR==1 & dat$female==0,])
t.test(Q14>=3 ~ Q14_hand1, data=dat[dat$Q14_HR==2 & dat$female==0,])

dat$Q14_hand1.b <- ifelse(dat$Q14_hand1==2, 1, 0)
summary(lm(Q14>=3 ~ Q14_hand1.b*female, data=dat[dat$Q14_HR==1,]))
summary(lm(Q14>=3 ~ Q14_hand1.b*female, data=dat[dat$Q14_HR==2,]))

summary(lm(Q14>=3 ~ Q14_hand1.b*dat$Q14_HR.b*female, data=dat))
dev.off()

##############################################################################################
### By IMCP
###############################################################################################
# Table A9: Logit Results: Should be Asked to Shake Hands [%]. Full Sample.

m.exp.1 <- glm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b*imcp, data=dat, family=binomial(link=logit))
sim.exp.1 <- sim(m.exp.1, 1000)

m.exp.2 <- glm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b*imcp, data=dat, family=binomial(link=logit))
sim.exp.2 <- sim(m.exp.2, 1000)

m.exp.3 <- glm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b*imcp, data=dat, family=binomial(link=logit))
sim.exp.3 <- sim(m.exp.3, 1000)

m.exp.1.b <- glm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b*imcp + female + age + abi, data=dat, family=binomial(link=logit))
m.exp.2.b <- glm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b*imcp + female + age + abi, data=dat, family=binomial(link=logit))
m.exp.3.b <- glm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b*imcp + female + age + abi, data=dat, family=binomial(link=logit))

stargazer(m.exp.1, m.exp.1.b, m.exp.2, m.exp.2.b, m.exp.3, m.exp.3.b, type="text", omit.stat=c("adj.rsq", "f", "ser", "ll", "aic"))

# Table A12: Logit Results: Should be Asked to Shake Hands [%]. Passed Attention Check.

m.exp.1 <- glm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b*imcp, data=dat2, family=binomial(link=logit))
m.exp.2 <- glm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b*imcp, data=dat2, family=binomial(link=logit))
m.exp.3 <- glm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b*imcp, data=dat2, family=binomial(link=logit))

m.exp.1.b <- glm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b*imcp + female + age + abi, data=dat2, family=binomial(link=logit))
m.exp.2.b <- glm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b*imcp + female + age + abi, data=dat2, family=binomial(link=logit))
m.exp.3.b <- glm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b*imcp + female + age + abi, data=dat2, family=binomial(link=logit))

stargazer(m.exp.1, m.exp.1.b, m.exp.2, m.exp.2.b, m.exp.3, m.exp.3.b, type="text", omit.stat=c("adj.rsq", "f", "ser", "ll", "aic"))

# Table A16: Brant test of null hypothesis that the ‘parallel regression’ assumption holds.

ord.1 <- polr(as.factor(Q12) ~ Q12_hand1.b*Q12_politician.b*imcp, data=dat2)
ord.2 <- polr(as.factor(Q13) ~ Q13_hand1.b*Q13_teacher.b*imcp, data=dat2)
ord.3 <- polr(as.factor(Q14) ~ Q14_hand1.b*Q14_HR.b*imcp, data=dat2)

brant(ord.1)
brant(ord.2)
brant(ord.3)

##############################################################################################
# Figure 3: Results of Handshaking Experiments by Respondents’ Internal Motivation to Control Prejudice (IMCP). 
png("bjps_handshaking_figure_3.png", width=6, height=9, units="in", res=600)

### Local Community Center
ci.m <- matrix(NA, 5, 4)
ci.f <- matrix(NA, 5, 4)
n <- 1
for(x in 1:5){
  
  # Male
  # Hand on Heart
  ci.m[n, 1] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2] + coef(sim.exp.1)[,3]*0 + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x +
                                    coef(sim.exp.1)[,7]*x*0 +  coef(sim.exp.1)[,8]*x*0), c(.025)) 
  
  ci.m[n, 2] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2] + coef(sim.exp.1)[,3]*0 + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x +
                                    coef(sim.exp.1)[,7]*x*0 +  coef(sim.exp.1)[,8]*x*0), c(.975)) 
  
  # No Mention
  ci.m[n, 3] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2]*0 + coef(sim.exp.1)[,3]*0 + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x*0 +
                                    coef(sim.exp.1)[,7]*x*0 +  coef(sim.exp.1)[,8]*x*0*0), c(.025)) 
  
  ci.m[n, 4] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2]*0 + coef(sim.exp.1)[,3]*0 + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x*0 +
                                    coef(sim.exp.1)[,7]*x*0 +  coef(sim.exp.1)[,8]*x*0*0), c(.975)) 
  
  # Female
  # Hand on Heart
  ci.f[n, 1] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2] + coef(sim.exp.1)[,3] + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5] + coef(sim.exp.1)[,6]*x +
                                    coef(sim.exp.1)[,7]*x +  coef(sim.exp.1)[,8]*x), c(.025)) 
  
  ci.f[n, 2] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2] + coef(sim.exp.1)[,3] + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5] + coef(sim.exp.1)[,6]*x +
                                    coef(sim.exp.1)[,7]*x +  coef(sim.exp.1)[,8]*x), c(.975)) 
  
  # No Mention
  ci.f[n, 3] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2]*0 + coef(sim.exp.1)[,3] + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x*0 +
                                    coef(sim.exp.1)[,7]*x +  coef(sim.exp.1)[,8]*x*0), c(.025)) 
  
  ci.f[n, 4] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2]*0 + coef(sim.exp.1)[,3] + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x*0 +
                                    coef(sim.exp.1)[,7]*x +  coef(sim.exp.1)[,8]*x*0), c(.975)) 
  
  n <- n + 1
}

par(mfrow=c(3,2), tck=-.02, font.main=1, mgp=c(2, .5, 0), mar=c(6,3,2,1), oma=c(1,1,3,1))
# All
curve(invlogit(coef(m.exp.1)[1] + coef(m.exp.1)[2] + coef(m.exp.1)[3]*0 + 
                 coef(m.exp.1)[4]*x + coef(m.exp.1)[5]*0 + coef(m.exp.1)[6]*x +
                 coef(m.exp.1)[7]*x*0 +  coef(m.exp.1)[8]*x*0), ylim=c(0, 1), xlim=c(1, 5), xlab="IMCP [1-5]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Politicians", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.1)[1] + coef(m.exp.1)[2]*0 + coef(m.exp.1)[3]*0 + 
                 coef(m.exp.1)[4]*x + coef(m.exp.1)[5]*0 + coef(m.exp.1)[6]*x*0 +
                 coef(m.exp.1)[7]*x*0 +  coef(m.exp.1)[8]*x*0*0), add=T, lty=2) # male

segments(1:5, ci.m[, 1], 1:5, ci.m[, 2])
segments(1:5, ci.m[, 3], 1:5, ci.m[, 4])

text(4, ci.m[4, 4], "No Mention", pos=3, cex=.8)
text(4, ci.m[4, 1], "Hand on Heart", pos=1, cex=.8)

# Female
curve(invlogit(coef(m.exp.1)[1] + coef(m.exp.1)[2] + coef(m.exp.1)[3] + 
                 coef(m.exp.1)[4]*x + coef(m.exp.1)[5] + coef(m.exp.1)[6]*x +
                 coef(m.exp.1)[7]*x +  coef(m.exp.1)[8]*x), ylim=c(0, 1), xlim=c(1, 5), xlab="IMCP [1-5]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Female Politicians", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.1)[1] + coef(m.exp.1)[2]*0 + coef(m.exp.1)[3] + 
                 coef(m.exp.1)[4]*x + coef(m.exp.1)[5]*0 + coef(m.exp.1)[6]*x*0 +
                 coef(m.exp.1)[7]*x +  coef(m.exp.1)[8]*x*0*0), add=T, lty=2) # male

segments(1:5, ci.f[, 1], 1:5, ci.f[, 2])
segments(1:5, ci.f[, 3], 1:5, ci.f[, 4])
mtext("Local Community Center", side=3, line=1, outer=T, adj=0)

text(4, ci.f[4, 4], "No Mention", pos=3, cex=.8)
text(4, ci.f[4, 1], "Hand on Heart", pos=1, cex=.8)


### Local School
ci.m <- matrix(NA, 5, 4)
ci.f <- matrix(NA, 5, 4)
n <- 1
for(x in 1:5){
  
  # Male
  # Hand on Heart
  ci.m[n, 1] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2] + coef(sim.exp.2)[,3]*0 + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x +
                                    coef(sim.exp.2)[,7]*x*0 +  coef(sim.exp.2)[,8]*x*0), c(.025)) 
  
  ci.m[n, 2] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2] + coef(sim.exp.2)[,3]*0 + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x +
                                    coef(sim.exp.2)[,7]*x*0 +  coef(sim.exp.2)[,8]*x*0), c(.975)) 
  
  # No Mention
  ci.m[n, 3] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2]*0 + coef(sim.exp.2)[,3]*0 + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x*0 +
                                    coef(sim.exp.2)[,7]*x*0 +  coef(sim.exp.2)[,8]*x*0*0), c(.025)) 
  
  ci.m[n, 4] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2]*0 + coef(sim.exp.2)[,3]*0 + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x*0 +
                                    coef(sim.exp.2)[,7]*x*0 +  coef(sim.exp.2)[,8]*x*0*0), c(.975)) 
  
  # Female
  # Hand on Heart
  ci.f[n, 1] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2] + coef(sim.exp.2)[,3] + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5] + coef(sim.exp.2)[,6]*x +
                                    coef(sim.exp.2)[,7]*x +  coef(sim.exp.2)[,8]*x), c(.025)) 
  
  ci.f[n, 2] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2] + coef(sim.exp.2)[,3] + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5] + coef(sim.exp.2)[,6]*x +
                                    coef(sim.exp.2)[,7]*x +  coef(sim.exp.2)[,8]*x), c(.975)) 
  
  # No Mention
  ci.f[n, 3] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2]*0 + coef(sim.exp.2)[,3] + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x*0 +
                                    coef(sim.exp.2)[,7]*x +  coef(sim.exp.2)[,8]*x*0), c(.025)) 
  
  ci.f[n, 4] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2]*0 + coef(sim.exp.2)[,3] + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x*0 +
                                    coef(sim.exp.2)[,7]*x +  coef(sim.exp.2)[,8]*x*0), c(.975)) 
  
  n <- n + 1
}

lines(seq(1,5, by=.1), ci.m[,1])
lines(seq(1,5, by=.1), ci.m[,2])
lines(seq(1,5, by=.1), ci.m[,3])
lines(seq(1,5, by=.1), ci.m[,4])



#par(mfrow=c(1,2), tck=-.02, font.main=1, mgp=c(2, .5, 0), mar=c(3,3,1,1), oma=c(1,1,2,1))
# All
curve(invlogit(coef(m.exp.2)[1] + coef(m.exp.2)[2] + coef(m.exp.2)[3]*0 + 
                 coef(m.exp.2)[4]*x + coef(m.exp.2)[5]*0 + coef(m.exp.2)[6]*x +
                 coef(m.exp.2)[7]*x*0 +  coef(m.exp.2)[8]*x*0), ylim=c(0, 1), xlim=c(1, 5), xlab="IMCP [1-5]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Teachers", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.2)[1] + coef(m.exp.2)[2]*0 + coef(m.exp.2)[3]*0 + 
                 coef(m.exp.2)[4]*x + coef(m.exp.2)[5]*0 + coef(m.exp.2)[6]*x*0 +
                 coef(m.exp.2)[7]*x*0 +  coef(m.exp.2)[8]*x*0*0), add=T, lty=2) # male

segments(1:5, ci.m[, 1], 1:5, ci.m[, 2])
segments(1:5, ci.m[, 3], 1:5, ci.m[, 4])

text(4, ci.m[4, 4], "No Mention", pos=3, cex=.8)
text(4, ci.m[4, 1], "Hand on Heart", pos=1, cex=.8)

# Female
curve(invlogit(coef(m.exp.2)[1] + coef(m.exp.2)[2] + coef(m.exp.2)[3] + 
                 coef(m.exp.2)[4]*x + coef(m.exp.2)[5] + coef(m.exp.2)[6]*x +
                 coef(m.exp.2)[7]*x +  coef(m.exp.2)[8]*x), ylim=c(0, 1), xlim=c(1, 5), xlab="IMCP [1-5]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Female Teachers", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.2)[1] + coef(m.exp.2)[2]*0 + coef(m.exp.2)[3] + 
                 coef(m.exp.2)[4]*x + coef(m.exp.2)[5]*0 + coef(m.exp.2)[6]*x*0 +
                 coef(m.exp.2)[7]*x +  coef(m.exp.2)[8]*x*0*0), add=T, lty=2) # male

segments(1:5, ci.f[, 1], 1:5, ci.f[, 2])
segments(1:5, ci.f[, 3], 1:5, ci.f[, 4])
mtext("Local School", side=3, line=-21, outer=T, adj=0)

text(4, ci.f[4, 4], "No Mention", pos=3, cex=.8)
text(4, ci.f[4, 1], "Hand on Heart", pos=1, cex=.8)

### Local Company
ci.m <- matrix(NA, 5, 4)
ci.f <- matrix(NA, 5, 4)
n <- 1
for(x in 1:5){
  
  # Male
  # Hand on Heart
  ci.m[n, 1] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2] + coef(sim.exp.3)[,3]*0 + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x +
                                    coef(sim.exp.3)[,7]*x*0 +  coef(sim.exp.3)[,8]*x*0), c(.025)) 
  
  ci.m[n, 2] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2] + coef(sim.exp.3)[,3]*0 + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x +
                                    coef(sim.exp.3)[,7]*x*0 +  coef(sim.exp.3)[,8]*x*0), c(.975)) 
  
  # No Mention
  ci.m[n, 3] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2]*0 + coef(sim.exp.3)[,3]*0 + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x*0 +
                                    coef(sim.exp.3)[,7]*x*0 +  coef(sim.exp.3)[,8]*x*0*0), c(.025)) 
  
  ci.m[n, 4] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2]*0 + coef(sim.exp.3)[,3]*0 + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x*0 +
                                    coef(sim.exp.3)[,7]*x*0 +  coef(sim.exp.3)[,8]*x*0*0), c(.975)) 
  
  # Female
  # Hand on Heart
  ci.f[n, 1] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2] + coef(sim.exp.3)[,3] + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5] + coef(sim.exp.3)[,6]*x +
                                    coef(sim.exp.3)[,7]*x +  coef(sim.exp.3)[,8]*x), c(.025)) 
  
  ci.f[n, 2] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2] + coef(sim.exp.3)[,3] + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5] + coef(sim.exp.3)[,6]*x +
                                    coef(sim.exp.3)[,7]*x +  coef(sim.exp.3)[,8]*x), c(.975)) 
  
  # No Mention
  ci.f[n, 3] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2]*0 + coef(sim.exp.3)[,3] + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x*0 +
                                    coef(sim.exp.3)[,7]*x +  coef(sim.exp.3)[,8]*x*0), c(.025)) 
  
  ci.f[n, 4] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2]*0 + coef(sim.exp.3)[,3] + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x*0 +
                                    coef(sim.exp.3)[,7]*x +  coef(sim.exp.3)[,8]*x*0), c(.975)) 
  
  n <- n + 1
}


#par(mfrow=c(1,2), tck=-.02, font.main=1, mgp=c(2, .5, 0), mar=c(3,3,1,1), oma=c(1,1,2,1))
# All
curve(invlogit(coef(m.exp.3)[1] + coef(m.exp.3)[2] + coef(m.exp.3)[3]*0 + 
                 coef(m.exp.3)[4]*x + coef(m.exp.3)[5]*0 + coef(m.exp.3)[6]*x +
                 coef(m.exp.3)[7]*x*0 +  coef(m.exp.3)[8]*x*0), ylim=c(0, 1), xlim=c(1, 5), xlab="IMCP [1-5]", 
      ylab="Should Be Asked to Shake Hands [%]", main="HR Managers", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.3)[1] + coef(m.exp.3)[2]*0 + coef(m.exp.3)[3]*0 + 
                 coef(m.exp.3)[4]*x + coef(m.exp.3)[5]*0 + coef(m.exp.3)[6]*x*0 +
                 coef(m.exp.3)[7]*x*0 +  coef(m.exp.3)[8]*x*0*0), add=T, lty=2) # male

segments(1:5, ci.m[, 1], 1:5, ci.m[, 2])
segments(1:5, ci.m[, 3], 1:5, ci.m[, 4])

text(4, ci.m[4, 4], "No Mention", pos=3, cex=.8)
text(4, ci.m[4, 1], "Hand on Heart", pos=1, cex=.8)

# Female
curve(invlogit(coef(m.exp.3)[1] + coef(m.exp.3)[2] + coef(m.exp.3)[3] + 
                 coef(m.exp.3)[4]*x + coef(m.exp.3)[5] + coef(m.exp.3)[6]*x +
                 coef(m.exp.3)[7]*x +  coef(m.exp.3)[8]*x), ylim=c(0, 1), xlim=c(1, 5), xlab="IMCP [1-5]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Female HR Managers", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.3)[1] + coef(m.exp.3)[2]*0 + coef(m.exp.3)[3] + 
                 coef(m.exp.3)[4]*x + coef(m.exp.3)[5]*0 + coef(m.exp.3)[6]*x*0 +
                 coef(m.exp.3)[7]*x +  coef(m.exp.3)[8]*x*0*0), add=T, lty=2) # male

segments(1:5, ci.f[, 1], 1:5, ci.f[, 2])
segments(1:5, ci.f[, 3], 1:5, ci.f[, 4])
mtext("Local Company", side=3, line=-42, outer=T, adj=0)

text(4, ci.f[4, 4], "No Mention", pos=3, cex=.8)
text(4, ci.f[4, 1], "Hand on Heart", pos=1, cex=.8)
dev.off()

###############################################################################################
### By Immigration Attitudes
###############################################################################################
# Table A10: Logit Results: Should be Asked to Shake Hands [%]. Full Sample.

m.exp.1 <- glm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b*immig, data=dat, family=binomial(link=logit))
sim.exp.1 <- sim(m.exp.1, 1000)

m.exp.2 <- glm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b*immig, data=dat, family=binomial(link=logit))
sim.exp.2 <- sim(m.exp.2, 1000)

m.exp.3 <- glm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b*immig, data=dat, family=binomial(link=logit))
sim.exp.3 <- sim(m.exp.3, 1000)

m.exp.1.b <- glm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b*immig + female + age + abi, data=dat, family=binomial(link=logit))
m.exp.2.b <- glm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b*immig + female + age + abi, data=dat, family=binomial(link=logit))
m.exp.3.b <- glm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b*immig + female + age + abi, data=dat, family=binomial(link=logit))

stargazer(m.exp.1, m.exp.1.b, m.exp.2, m.exp.2.b, m.exp.3, m.exp.3.b, type="text", omit.stat=c("adj.rsq", "f", "ser", "ll", "aic"))

# Table A13: Logit Results: Should be Asked to Shake Hands [%]. Passed Attention Check.

m.exp.1 <- glm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b*immig, data=dat2, family=binomial(link=logit))
m.exp.2 <- glm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b*immig, data=dat2, family=binomial(link=logit))
m.exp.3 <- glm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b*immig, data=dat2, family=binomial(link=logit))

m.exp.1.b <- glm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b*immig + female + age + abi, data=dat2, family=binomial(link=logit))
m.exp.2.b <- glm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b*immig + female + age + abi, data=dat2, family=binomial(link=logit))
m.exp.3.b <- glm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b*immig + female + age + abi, data=dat2, family=binomial(link=logit))

stargazer(m.exp.1, m.exp.1.b, m.exp.2, m.exp.2.b, m.exp.3, m.exp.3.b, type="text", omit.stat=c("adj.rsq", "f", "ser", "ll", "aic"))

# Table A17: Brant test of null hypothesis that the ‘parallel regression’ assumption holds.

ord.1 <- polr(as.factor(Q12) ~ Q12_hand1.b*Q12_politician.b*immig, data=dat2)
ord.2 <- polr(as.factor(Q13) ~ Q13_hand1.b*Q13_teacher.b*immig, data=dat2)
ord.3 <- polr(as.factor(Q14) ~ Q14_hand1.b*Q14_HR.b*immig, data=dat2)

brant(ord.1)
brant(ord.2)
brant(ord.3)

###############################################################################################
# Figure 4: Results of Handshaking Experiments by Respondents’ Immigration Attitude.  
png("bjps_handshaking_figure_4.png", width=6, height=9, units="in", res=600)

### Local Community Center
ci.m <- matrix(NA, 11, 4)
ci.f <- matrix(NA, 11, 4)
n <- 1
for(x in 0:10){
  
  # Male
  # Hand on Heart
  ci.m[n, 1] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2] + coef(sim.exp.1)[,3]*0 + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x +
                                    coef(sim.exp.1)[,7]*x*0 +  coef(sim.exp.1)[,8]*x*0), c(.025)) 
  
  ci.m[n, 2] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2] + coef(sim.exp.1)[,3]*0 + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x +
                                    coef(sim.exp.1)[,7]*x*0 +  coef(sim.exp.1)[,8]*x*0), c(.975)) 
  
  # No Mention
  ci.m[n, 3] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2]*0 + coef(sim.exp.1)[,3]*0 + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x*0 +
                                    coef(sim.exp.1)[,7]*x*0 +  coef(sim.exp.1)[,8]*x*0*0), c(.025)) 
  
  ci.m[n, 4] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2]*0 + coef(sim.exp.1)[,3]*0 + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x*0 +
                                    coef(sim.exp.1)[,7]*x*0 +  coef(sim.exp.1)[,8]*x*0*0), c(.975)) 
  
  # Female
  # Hand on Heart
  ci.f[n, 1] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2] + coef(sim.exp.1)[,3] + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5] + coef(sim.exp.1)[,6]*x +
                                    coef(sim.exp.1)[,7]*x +  coef(sim.exp.1)[,8]*x), c(.025)) 
  
  ci.f[n, 2] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2] + coef(sim.exp.1)[,3] + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5] + coef(sim.exp.1)[,6]*x +
                                    coef(sim.exp.1)[,7]*x +  coef(sim.exp.1)[,8]*x), c(.975)) 
  
  # No Mention
  ci.f[n, 3] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2]*0 + coef(sim.exp.1)[,3] + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x*0 +
                                    coef(sim.exp.1)[,7]*x +  coef(sim.exp.1)[,8]*x*0), c(.025)) 
  
  ci.f[n, 4] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2]*0 + coef(sim.exp.1)[,3] + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x*0 +
                                    coef(sim.exp.1)[,7]*x +  coef(sim.exp.1)[,8]*x*0), c(.975)) 
  
  n <- n + 1
}


par(mfrow=c(3,2), tck=-.02, font.main=1, mgp=c(2, .5, 0), mar=c(6,3,1,1), oma=c(1,1,3,1))
# All
curve(invlogit(coef(m.exp.1)[1] + coef(m.exp.1)[2] + coef(m.exp.1)[3]*0 + 
                 coef(m.exp.1)[4]*x + coef(m.exp.1)[5]*0 + coef(m.exp.1)[6]*x +
                 coef(m.exp.1)[7]*x*0 +  coef(m.exp.1)[8]*x*0), ylim=c(0, 1), xlim=c(0, 10), xlab="Immigration Attitude [0-10]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Politicians", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.1)[1] + coef(m.exp.1)[2]*0 + coef(m.exp.1)[3]*0 + 
                 coef(m.exp.1)[4]*x + coef(m.exp.1)[5]*0 + coef(m.exp.1)[6]*x*0 +
                 coef(m.exp.1)[7]*x*0 +  coef(m.exp.1)[8]*x*0*0), add=T, lty=2) # male

segments(0:10, ci.m[, 1], 0:10, ci.m[, 2])
segments(0:10, ci.m[, 3], 0:10, ci.m[, 4])

text(8, ci.m[9, 4], "No Mention", pos=3, cex=.8)
text(8, ci.m[9, 1], "Hand on Heart", pos=1, cex=.8)

# Female
curve(invlogit(coef(m.exp.1)[1] + coef(m.exp.1)[2] + coef(m.exp.1)[3] + 
                 coef(m.exp.1)[4]*x + coef(m.exp.1)[5] + coef(m.exp.1)[6]*x +
                 coef(m.exp.1)[7]*x +  coef(m.exp.1)[8]*x), ylim=c(0, 1), xlim=c(0, 10), xlab="Immigration Attitude [0-10]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Female Politicians", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.1)[1] + coef(m.exp.1)[2]*0 + coef(m.exp.1)[3] + 
                 coef(m.exp.1)[4]*x + coef(m.exp.1)[5]*0 + coef(m.exp.1)[6]*x*0 +
                 coef(m.exp.1)[7]*x +  coef(m.exp.1)[8]*x*0*0), add=T, lty=2) # male

segments(0:10, ci.f[, 1], 0:10, ci.f[, 2])
segments(0:10, ci.f[, 3], 0:10, ci.f[, 4])
mtext("Local Community Center", side=3, line=1, outer=T, adj=0)

text(8, ci.f[9, 4], "No Mention", pos=3, cex=.8)
text(8, ci.f[9, 1], "Hand on Heart", pos=1, cex=.8)

### Local School
ci.m <- matrix(NA, 11, 4)
ci.f <- matrix(NA, 11, 4)
n <- 1
for(x in 0:11){
  
  # Male
  # Hand on Heart
  ci.m[n, 1] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2] + coef(sim.exp.2)[,3]*0 + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x +
                                    coef(sim.exp.2)[,7]*x*0 +  coef(sim.exp.2)[,8]*x*0), c(.025)) 
  
  ci.m[n, 2] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2] + coef(sim.exp.2)[,3]*0 + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x +
                                    coef(sim.exp.2)[,7]*x*0 +  coef(sim.exp.2)[,8]*x*0), c(.975)) 
  
  # No Mention
  ci.m[n, 3] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2]*0 + coef(sim.exp.2)[,3]*0 + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x*0 +
                                    coef(sim.exp.2)[,7]*x*0 +  coef(sim.exp.2)[,8]*x*0*0), c(.025)) 
  
  ci.m[n, 4] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2]*0 + coef(sim.exp.2)[,3]*0 + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x*0 +
                                    coef(sim.exp.2)[,7]*x*0 +  coef(sim.exp.2)[,8]*x*0*0), c(.975)) 
  
  # Female
  # Hand on Heart
  ci.f[n, 1] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2] + coef(sim.exp.2)[,3] + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5] + coef(sim.exp.2)[,6]*x +
                                    coef(sim.exp.2)[,7]*x +  coef(sim.exp.2)[,8]*x), c(.025)) 
  
  ci.f[n, 2] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2] + coef(sim.exp.2)[,3] + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5] + coef(sim.exp.2)[,6]*x +
                                    coef(sim.exp.2)[,7]*x +  coef(sim.exp.2)[,8]*x), c(.975)) 
  
  # No Mention
  ci.f[n, 3] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2]*0 + coef(sim.exp.2)[,3] + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x*0 +
                                    coef(sim.exp.2)[,7]*x +  coef(sim.exp.2)[,8]*x*0), c(.025)) 
  
  ci.f[n, 4] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2]*0 + coef(sim.exp.2)[,3] + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x*0 +
                                    coef(sim.exp.2)[,7]*x +  coef(sim.exp.2)[,8]*x*0), c(.975)) 
  
  n <- n + 1
}


#par(mfrow=c(1,2), tck=-.02, font.main=1, mgp=c(2, .5, 0), mar=c(3,3,1,1), oma=c(1,1,2,1))
# All
curve(invlogit(coef(m.exp.2)[1] + coef(m.exp.2)[2] + coef(m.exp.2)[3]*0 + 
                 coef(m.exp.2)[4]*x + coef(m.exp.2)[5]*0 + coef(m.exp.2)[6]*x +
                 coef(m.exp.2)[7]*x*0 +  coef(m.exp.2)[8]*x*0), ylim=c(0, 1), xlim=c(0, 10), xlab="Immigration Attitude [0-10]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Teachers", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.2)[1] + coef(m.exp.2)[2]*0 + coef(m.exp.2)[3]*0 + 
                 coef(m.exp.2)[4]*x + coef(m.exp.2)[5]*0 + coef(m.exp.2)[6]*x*0 +
                 coef(m.exp.2)[7]*x*0 +  coef(m.exp.2)[8]*x*0*0), add=T, lty=2) # male

segments(0:10, ci.m[, 1], 0:10, ci.m[, 2])
segments(0:10, ci.m[, 3], 0:10, ci.m[, 4])

text(8, ci.m[9, 4], "No Mention", pos=3, cex=.8)
text(8, ci.m[9, 1], "Hand on Heart", pos=1, cex=.8)

# Female
curve(invlogit(coef(m.exp.2)[1] + coef(m.exp.2)[2] + coef(m.exp.2)[3] + 
                 coef(m.exp.2)[4]*x + coef(m.exp.2)[5] + coef(m.exp.2)[6]*x +
                 coef(m.exp.2)[7]*x +  coef(m.exp.2)[8]*x), ylim=c(0, 1), xlim=c(0, 10), xlab="Immigration Attitude [0-10]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Female Teachers", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.2)[1] + coef(m.exp.2)[2]*0 + coef(m.exp.2)[3] + 
                 coef(m.exp.2)[4]*x + coef(m.exp.2)[5]*0 + coef(m.exp.2)[6]*x*0 +
                 coef(m.exp.2)[7]*x +  coef(m.exp.2)[8]*x*0*0), add=T, lty=2) # male

segments(0:10, ci.f[, 1], 0:10, ci.f[, 2])
segments(0:10, ci.f[, 3], 0:10, ci.f[, 4])
mtext("Local School", side=3, line=-21, outer=T, adj=0)

text(8, ci.f[9, 4], "No Mention", pos=3, cex=.8)
text(8, ci.f[9, 1], "Hand on Heart", pos=1, cex=.8)


### Local Company
ci.m <- matrix(NA, 11, 4)
ci.f <- matrix(NA, 11, 4)
n <- 1
for(x in 0:10){
  
  # Male
  # Hand on Heart
  ci.m[n, 1] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2] + coef(sim.exp.3)[,3]*0 + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x +
                                    coef(sim.exp.3)[,7]*x*0 +  coef(sim.exp.3)[,8]*x*0), c(.025)) 
  
  ci.m[n, 2] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2] + coef(sim.exp.3)[,3]*0 + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x +
                                    coef(sim.exp.3)[,7]*x*0 +  coef(sim.exp.3)[,8]*x*0), c(.975)) 
  
  # No Mention
  ci.m[n, 3] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2]*0 + coef(sim.exp.3)[,3]*0 + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x*0 +
                                    coef(sim.exp.3)[,7]*x*0 +  coef(sim.exp.3)[,8]*x*0*0), c(.025)) 
  
  ci.m[n, 4] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2]*0 + coef(sim.exp.3)[,3]*0 + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x*0 +
                                    coef(sim.exp.3)[,7]*x*0 +  coef(sim.exp.3)[,8]*x*0*0), c(.975)) 
  
  # Female
  # Hand on Heart
  ci.f[n, 1] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2] + coef(sim.exp.3)[,3] + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5] + coef(sim.exp.3)[,6]*x +
                                    coef(sim.exp.3)[,7]*x +  coef(sim.exp.3)[,8]*x), c(.025)) 
  
  ci.f[n, 2] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2] + coef(sim.exp.3)[,3] + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5] + coef(sim.exp.3)[,6]*x +
                                    coef(sim.exp.3)[,7]*x +  coef(sim.exp.3)[,8]*x), c(.975)) 
  
  # No Mention
  ci.f[n, 3] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2]*0 + coef(sim.exp.3)[,3] + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x*0 +
                                    coef(sim.exp.3)[,7]*x +  coef(sim.exp.3)[,8]*x*0), c(.025)) 
  
  ci.f[n, 4] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2]*0 + coef(sim.exp.3)[,3] + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x*0 +
                                    coef(sim.exp.3)[,7]*x +  coef(sim.exp.3)[,8]*x*0), c(.975)) 
  
  n <- n + 1
}

#par(mfrow=c(1,2), tck=-.02, font.main=1, mgp=c(2, .5, 0), mar=c(3,3,1,1), oma=c(1,1,2,1))
# All
curve(invlogit(coef(m.exp.3)[1] + coef(m.exp.3)[2] + coef(m.exp.3)[3]*0 + 
                 coef(m.exp.3)[4]*x + coef(m.exp.3)[5]*0 + coef(m.exp.3)[6]*x +
                 coef(m.exp.3)[7]*x*0 +  coef(m.exp.3)[8]*x*0), ylim=c(0, 1), xlim=c(0, 10), xlab="Immigration Attitude [0-10]", 
      ylab="Should Be Asked to Shake Hands [%]", main="HR Managers", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.3)[1] + coef(m.exp.3)[2]*0 + coef(m.exp.3)[3]*0 + 
                 coef(m.exp.3)[4]*x + coef(m.exp.3)[5]*0 + coef(m.exp.3)[6]*x*0 +
                 coef(m.exp.3)[7]*x*0 +  coef(m.exp.3)[8]*x*0*0), add=T, lty=2) # male

segments(0:10, ci.m[, 1], 0:10, ci.m[, 2])
segments(0:10, ci.m[, 3], 0:10, ci.m[, 4])

text(8, ci.m[9, 4], "No Mention", pos=3, cex=.8)
text(8, ci.m[9, 1], "Hand on Heart", pos=1, cex=.8)

# Female
curve(invlogit(coef(m.exp.3)[1] + coef(m.exp.3)[2] + coef(m.exp.3)[3] + 
                 coef(m.exp.3)[4]*x + coef(m.exp.3)[5] + coef(m.exp.3)[6]*x +
                 coef(m.exp.3)[7]*x +  coef(m.exp.3)[8]*x), ylim=c(0, 1), xlim=c(0, 10), xlab="Immigration Attitude [0-10]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Female HR Managers", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.3)[1] + coef(m.exp.3)[2]*0 + coef(m.exp.3)[3] + 
                 coef(m.exp.3)[4]*x + coef(m.exp.3)[5]*0 + coef(m.exp.3)[6]*x*0 +
                 coef(m.exp.3)[7]*x +  coef(m.exp.3)[8]*x*0*0), add=T, lty=2) # male

segments(0:10, ci.f[, 1], 0:10, ci.f[, 2])
segments(0:10, ci.f[, 3], 0:10, ci.f[, 4])
mtext("Local Company", side=3, line=-42, outer=T, adj=0)

text(8, ci.f[9, 4], "No Mention", pos=3, cex=.8)
text(8, ci.f[9, 1], "Hand on Heart", pos=1, cex=.8)
dev.off()

###############################################################################################
### By Left-Right-Ideology
###############################################################################################
# Table A11: Logit Results: Should be Asked to Shake Hands [%]. Full Sample.

m.exp.1 <- glm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b*leftright, data=dat, family=binomial(link=logit))
sim.exp.1 <- sim(m.exp.1, 1000)

m.exp.2 <- glm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b*leftright, data=dat, family=binomial(link=logit))
sim.exp.2 <- sim(m.exp.2, 1000)

m.exp.3 <- glm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b*leftright, data=dat, family=binomial(link=logit))
sim.exp.3 <- sim(m.exp.3, 1000)

m.exp.1.b <- glm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b*leftright + female + age + abi, data=dat, family=binomial(link=logit))
m.exp.2.b <- glm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b*leftright + female + age + abi, data=dat, family=binomial(link=logit))
m.exp.3.b <- glm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b*leftright + female + age + abi, data=dat, family=binomial(link=logit))

stargazer(m.exp.1, m.exp.1.b, m.exp.2, m.exp.2.b, m.exp.3, m.exp.3.b, type="text", omit.stat=c("adj.rsq", "f", "ser", "ll", "aic"))

# Table A14: Logit Results: Should be Asked to Shake Hands [%]. Passed the Attention Check

m.exp.1 <- glm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b*leftright, data=dat2, family=binomial(link=logit))
m.exp.2 <- glm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b*leftright, data=dat2, family=binomial(link=logit))
m.exp.3 <- glm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b*leftright, data=dat2, family=binomial(link=logit))

m.exp.1.b <- glm(Q12>=3 ~ Q12_hand1.b*Q12_politician.b*leftright + female + age + abi, data=dat2, family=binomial(link=logit))
m.exp.2.b <- glm(Q13>=3 ~ Q13_hand1.b*Q13_teacher.b*leftright + female + age + abi, data=dat2, family=binomial(link=logit))
m.exp.3.b <- glm(Q14>=3 ~ Q14_hand1.b*Q14_HR.b*leftright + female + age + abi, data=dat2, family=binomial(link=logit))

stargazer(m.exp.1, m.exp.1.b, m.exp.2, m.exp.2.b, m.exp.3, m.exp.3.b, type="text", omit.stat=c("adj.rsq", "f", "ser", "ll", "aic"))

# Table A18: Brant test of null hypothesis that the ‘parallel regression’ assumption holds.

ord.1 <- polr(as.factor(Q12) ~ Q12_hand1.b*Q12_politician.b*leftright, data=dat2)
ord.2 <- polr(as.factor(Q13) ~ Q13_hand1.b*Q13_teacher.b*leftright, data=dat2)
ord.3 <- polr(as.factor(Q14) ~ Q14_hand1.b*Q14_HR.b*leftright, data=dat2)

brant(ord.1)
brant(ord.2)
brant(ord.3)


###############################################################################################
# Figure 5: Results of Handshaking Experiments by Respondents’ Political Left-Right-Ideology. 
png("bjps_handshaking_figure_5.png", width=6, height=9, units="in", res=600)

### Local Community Center
ci.m <- matrix(NA, 11, 4)
ci.f <- matrix(NA, 11, 4)
n <- 1
for(x in 0:10){
  
  # Male
  # Hand on Heart
  ci.m[n, 1] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2] + coef(sim.exp.1)[,3]*0 + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x +
                                    coef(sim.exp.1)[,7]*x*0 +  coef(sim.exp.1)[,8]*x*0), c(.025)) 
  
  ci.m[n, 2] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2] + coef(sim.exp.1)[,3]*0 + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x +
                                    coef(sim.exp.1)[,7]*x*0 +  coef(sim.exp.1)[,8]*x*0), c(.975)) 
  
  # No Mention
  ci.m[n, 3] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2]*0 + coef(sim.exp.1)[,3]*0 + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x*0 +
                                    coef(sim.exp.1)[,7]*x*0 +  coef(sim.exp.1)[,8]*x*0*0), c(.025)) 
  
  ci.m[n, 4] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2]*0 + coef(sim.exp.1)[,3]*0 + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x*0 +
                                    coef(sim.exp.1)[,7]*x*0 +  coef(sim.exp.1)[,8]*x*0*0), c(.975)) 
  
  # Female
  # Hand on Heart
  ci.f[n, 1] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2] + coef(sim.exp.1)[,3] + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5] + coef(sim.exp.1)[,6]*x +
                                    coef(sim.exp.1)[,7]*x +  coef(sim.exp.1)[,8]*x), c(.025)) 
  
  ci.f[n, 2] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2] + coef(sim.exp.1)[,3] + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5] + coef(sim.exp.1)[,6]*x +
                                    coef(sim.exp.1)[,7]*x +  coef(sim.exp.1)[,8]*x), c(.975)) 
  
  # No Mention
  ci.f[n, 3] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2]*0 + coef(sim.exp.1)[,3] + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x*0 +
                                    coef(sim.exp.1)[,7]*x +  coef(sim.exp.1)[,8]*x*0), c(.025)) 
  
  ci.f[n, 4] <- quantile(invlogit(coef(sim.exp.1)[,1] + coef(sim.exp.1)[,2]*0 + coef(sim.exp.1)[,3] + 
                                    coef(sim.exp.1)[,4]*x + coef(sim.exp.1)[,5]*0 + coef(sim.exp.1)[,6]*x*0 +
                                    coef(sim.exp.1)[,7]*x +  coef(sim.exp.1)[,8]*x*0), c(.975)) 
  
  n <- n + 1
}


par(mfrow=c(3,2), tck=-.02, font.main=1, mgp=c(2, .5, 0), mar=c(6,3,1,1), oma=c(1,1,3,1))
# All
curve(invlogit(coef(m.exp.1)[1] + coef(m.exp.1)[2] + coef(m.exp.1)[3]*0 + 
                 coef(m.exp.1)[4]*x + coef(m.exp.1)[5]*0 + coef(m.exp.1)[6]*x +
                 coef(m.exp.1)[7]*x*0 +  coef(m.exp.1)[8]*x*0), ylim=c(0, 1), xlim=c(0, 10), xlab="Left-Right-Ideology [0-10]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Politicians", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.1)[1] + coef(m.exp.1)[2]*0 + coef(m.exp.1)[3]*0 + 
                 coef(m.exp.1)[4]*x + coef(m.exp.1)[5]*0 + coef(m.exp.1)[6]*x*0 +
                 coef(m.exp.1)[7]*x*0 +  coef(m.exp.1)[8]*x*0*0), add=T, lty=2) # male

segments(0:10, ci.m[, 1], 0:10, ci.m[, 2])
segments(0:10, ci.m[, 3], 0:10, ci.m[, 4])

text(4, ci.m[4, 4], "No Mention", pos=3, cex=.8)
text(4, ci.m[4, 1], "Hand on Heart", pos=1, cex=.8)

# Female
curve(invlogit(coef(m.exp.1)[1] + coef(m.exp.1)[2] + coef(m.exp.1)[3] + 
                 coef(m.exp.1)[4]*x + coef(m.exp.1)[5] + coef(m.exp.1)[6]*x +
                 coef(m.exp.1)[7]*x +  coef(m.exp.1)[8]*x), ylim=c(0, 1), xlim=c(0, 10), xlab="Left-Right-Ideology [0-10]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Female Politicians", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.1)[1] + coef(m.exp.1)[2]*0 + coef(m.exp.1)[3] + 
                 coef(m.exp.1)[4]*x + coef(m.exp.1)[5]*0 + coef(m.exp.1)[6]*x*0 +
                 coef(m.exp.1)[7]*x +  coef(m.exp.1)[8]*x*0*0), add=T, lty=2) # male

segments(0:10, ci.f[, 1], 0:10, ci.f[, 2])
segments(0:10, ci.f[, 3], 0:10, ci.f[, 4])
mtext("Local Community Center", side=3, line=1, outer=T, adj=0)

text(4, ci.f[4, 4], "No Mention", pos=3, cex=.8)
text(4, ci.f[4, 1], "Hand on Heart", pos=1, cex=.8)

### Local School
ci.m <- matrix(NA, 11, 4)
ci.f <- matrix(NA, 11, 4)
n <- 1
for(x in 0:11){
  
  # Male
  # Hand on Heart
  ci.m[n, 1] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2] + coef(sim.exp.2)[,3]*0 + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x +
                                    coef(sim.exp.2)[,7]*x*0 +  coef(sim.exp.2)[,8]*x*0), c(.025)) 
  
  ci.m[n, 2] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2] + coef(sim.exp.2)[,3]*0 + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x +
                                    coef(sim.exp.2)[,7]*x*0 +  coef(sim.exp.2)[,8]*x*0), c(.975)) 
  
  # No Mention
  ci.m[n, 3] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2]*0 + coef(sim.exp.2)[,3]*0 + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x*0 +
                                    coef(sim.exp.2)[,7]*x*0 +  coef(sim.exp.2)[,8]*x*0*0), c(.025)) 
  
  ci.m[n, 4] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2]*0 + coef(sim.exp.2)[,3]*0 + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x*0 +
                                    coef(sim.exp.2)[,7]*x*0 +  coef(sim.exp.2)[,8]*x*0*0), c(.975)) 
  
  # Female
  # Hand on Heart
  ci.f[n, 1] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2] + coef(sim.exp.2)[,3] + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5] + coef(sim.exp.2)[,6]*x +
                                    coef(sim.exp.2)[,7]*x +  coef(sim.exp.2)[,8]*x), c(.025)) 
  
  ci.f[n, 2] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2] + coef(sim.exp.2)[,3] + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5] + coef(sim.exp.2)[,6]*x +
                                    coef(sim.exp.2)[,7]*x +  coef(sim.exp.2)[,8]*x), c(.975)) 
  
  # No Mention
  ci.f[n, 3] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2]*0 + coef(sim.exp.2)[,3] + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x*0 +
                                    coef(sim.exp.2)[,7]*x +  coef(sim.exp.2)[,8]*x*0), c(.025)) 
  
  ci.f[n, 4] <- quantile(invlogit(coef(sim.exp.2)[,1] + coef(sim.exp.2)[,2]*0 + coef(sim.exp.2)[,3] + 
                                    coef(sim.exp.2)[,4]*x + coef(sim.exp.2)[,5]*0 + coef(sim.exp.2)[,6]*x*0 +
                                    coef(sim.exp.2)[,7]*x +  coef(sim.exp.2)[,8]*x*0), c(.975)) 
  
  n <- n + 1
}


#par(mfrow=c(1,2), tck=-.02, font.main=1, mgp=c(2, .5, 0), mar=c(3,3,1,1), oma=c(1,1,2,1))
# All
curve(invlogit(coef(m.exp.2)[1] + coef(m.exp.2)[2] + coef(m.exp.2)[3]*0 + 
                 coef(m.exp.2)[4]*x + coef(m.exp.2)[5]*0 + coef(m.exp.2)[6]*x +
                 coef(m.exp.2)[7]*x*0 +  coef(m.exp.2)[8]*x*0), ylim=c(0, 1), xlim=c(0, 10), xlab="Left-Right-Ideology [0-10]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Teachers", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.2)[1] + coef(m.exp.2)[2]*0 + coef(m.exp.2)[3]*0 + 
                 coef(m.exp.2)[4]*x + coef(m.exp.2)[5]*0 + coef(m.exp.2)[6]*x*0 +
                 coef(m.exp.2)[7]*x*0 +  coef(m.exp.2)[8]*x*0*0), add=T, lty=2) # male

segments(0:10, ci.m[, 1], 0:10, ci.m[, 2])
segments(0:10, ci.m[, 3], 0:10, ci.m[, 4])

text(4, ci.m[4, 4], "No Mention", pos=3, cex=.8)
text(4, ci.m[4, 1], "Hand on Heart", pos=1, cex=.8)

# Female
curve(invlogit(coef(m.exp.2)[1] + coef(m.exp.2)[2] + coef(m.exp.2)[3] + 
                 coef(m.exp.2)[4]*x + coef(m.exp.2)[5] + coef(m.exp.2)[6]*x +
                 coef(m.exp.2)[7]*x +  coef(m.exp.2)[8]*x), ylim=c(0, 1), xlim=c(0, 10), xlab="Left-Right-Ideology [0-10]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Female Teachers", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.2)[1] + coef(m.exp.2)[2]*0 + coef(m.exp.2)[3] + 
                 coef(m.exp.2)[4]*x + coef(m.exp.2)[5]*0 + coef(m.exp.2)[6]*x*0 +
                 coef(m.exp.2)[7]*x +  coef(m.exp.2)[8]*x*0*0), add=T, lty=2) # male

segments(0:10, ci.f[, 1], 0:10, ci.f[, 2])
segments(0:10, ci.f[, 3], 0:10, ci.f[, 4])
mtext("Local School", side=3, line=-21, outer=T, adj=0)

text(4, ci.f[4, 4], "No Mention", pos=3, cex=.8)
text(4, ci.f[4, 1], "Hand on Heart", pos=1, cex=.8)


### Local Company
ci.m <- matrix(NA, 11, 4)
ci.f <- matrix(NA, 11, 4)
n <- 1
for(x in 0:10){
  
  # Male
  # Hand on Heart
  ci.m[n, 1] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2] + coef(sim.exp.3)[,3]*0 + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x +
                                    coef(sim.exp.3)[,7]*x*0 +  coef(sim.exp.3)[,8]*x*0), c(.025)) 
  
  ci.m[n, 2] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2] + coef(sim.exp.3)[,3]*0 + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x +
                                    coef(sim.exp.3)[,7]*x*0 +  coef(sim.exp.3)[,8]*x*0), c(.975)) 
  
  # No Mention
  ci.m[n, 3] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2]*0 + coef(sim.exp.3)[,3]*0 + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x*0 +
                                    coef(sim.exp.3)[,7]*x*0 +  coef(sim.exp.3)[,8]*x*0*0), c(.025)) 
  
  ci.m[n, 4] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2]*0 + coef(sim.exp.3)[,3]*0 + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x*0 +
                                    coef(sim.exp.3)[,7]*x*0 +  coef(sim.exp.3)[,8]*x*0*0), c(.975)) 
  
  # Female
  # Hand on Heart
  ci.f[n, 1] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2] + coef(sim.exp.3)[,3] + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5] + coef(sim.exp.3)[,6]*x +
                                    coef(sim.exp.3)[,7]*x +  coef(sim.exp.3)[,8]*x), c(.025)) 
  
  ci.f[n, 2] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2] + coef(sim.exp.3)[,3] + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5] + coef(sim.exp.3)[,6]*x +
                                    coef(sim.exp.3)[,7]*x +  coef(sim.exp.3)[,8]*x), c(.975)) 
  
  # No Mention
  ci.f[n, 3] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2]*0 + coef(sim.exp.3)[,3] + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x*0 +
                                    coef(sim.exp.3)[,7]*x +  coef(sim.exp.3)[,8]*x*0), c(.025)) 
  
  ci.f[n, 4] <- quantile(invlogit(coef(sim.exp.3)[,1] + coef(sim.exp.3)[,2]*0 + coef(sim.exp.3)[,3] + 
                                    coef(sim.exp.3)[,4]*x + coef(sim.exp.3)[,5]*0 + coef(sim.exp.3)[,6]*x*0 +
                                    coef(sim.exp.3)[,7]*x +  coef(sim.exp.3)[,8]*x*0), c(.975)) 
  
  n <- n + 1
}

#par(mfrow=c(1,2), tck=-.02, font.main=1, mgp=c(2, .5, 0), mar=c(3,3,1,1), oma=c(1,1,2,1))
# All
curve(invlogit(coef(m.exp.3)[1] + coef(m.exp.3)[2] + coef(m.exp.3)[3]*0 + 
                 coef(m.exp.3)[4]*x + coef(m.exp.3)[5]*0 + coef(m.exp.3)[6]*x +
                 coef(m.exp.3)[7]*x*0 +  coef(m.exp.3)[8]*x*0), ylim=c(0, 1), xlim=c(0, 10), xlab="Left-Right-Ideology [0-10]", 
      ylab="Should Be Asked to Shake Hands [%]", main="HR Managers", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.3)[1] + coef(m.exp.3)[2]*0 + coef(m.exp.3)[3]*0 + 
                 coef(m.exp.3)[4]*x + coef(m.exp.3)[5]*0 + coef(m.exp.3)[6]*x*0 +
                 coef(m.exp.3)[7]*x*0 +  coef(m.exp.3)[8]*x*0*0), add=T, lty=2) # male

segments(0:10, ci.m[, 1], 0:10, ci.m[, 2])
segments(0:10, ci.m[, 3], 0:10, ci.m[, 4])

text(4, ci.m[4, 4], "No Mention", pos=3, cex=.8)
text(4, ci.m[4, 1], "Hand on Heart", pos=1, cex=.8)

# Female
curve(invlogit(coef(m.exp.3)[1] + coef(m.exp.3)[2] + coef(m.exp.3)[3] + 
                 coef(m.exp.3)[4]*x + coef(m.exp.3)[5] + coef(m.exp.3)[6]*x +
                 coef(m.exp.3)[7]*x +  coef(m.exp.3)[8]*x), ylim=c(0, 1), xlim=c(0, 10), xlab="Left-Right-Ideology [0-10]", 
      ylab="Should Be Asked to Shake Hands [%]", main="Female HR Managers", axes=F) # male
abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))
axis(1)
box()

curve(invlogit(coef(m.exp.3)[1] + coef(m.exp.3)[2]*0 + coef(m.exp.3)[3] + 
                 coef(m.exp.3)[4]*x + coef(m.exp.3)[5]*0 + coef(m.exp.3)[6]*x*0 +
                 coef(m.exp.3)[7]*x +  coef(m.exp.3)[8]*x*0*0), add=T, lty=2) # male

segments(0:10, ci.f[, 1], 0:10, ci.f[, 2])
segments(0:10, ci.f[, 3], 0:10, ci.f[, 4])
mtext("Local Company", side=3, line=-42, outer=T, adj=0)

text(4, ci.f[4, 4], "No Mention", pos=3, cex=.8)
text(4, ci.f[4, 1], "Hand on Heart", pos=1, cex=.8)
dev.off()

###############################################################################################
# Replication: Symmetric Gender

load("bjps_handshaking_replication_data_gender_rep.RData")

dat$handshake_experiment <- ifelse(dat$handshake_experiment=="-99", NA, dat$handshake_experiment)
dat$handshake_experiment <- ifelse(dat$handshake_experiment=="rather_in_favor" | dat$handshake_experiment=="strongly_in_favor", 1, 0)

# Table A15: Replication of Handshaking Experiment Where Muslim Women Refuse to Shake Hands with Men: Should be Asked to Shake Hands [%]
means <- aggregate(dat$handshake_experiment, by=list(dat$handshake1, dat$handshake2), mean, na.rm=T)
sds <- aggregate(dat$handshake_experiment, by=list(dat$handshake1, dat$handshake2), sd, na.rm=T)
means
sds
table(dat$handshake1[is.na(dat$handshake_experiment)==F], dat$handshake2[is.na(dat$handshake_experiment)==F])

m <- lm(handshake_experiment ~ handshake1*handshake2, data=dat)

#############################################################################################
# Replication: Norway

load("bjps_handshaking_replication_data_norway_rep.RData")

norw$handshake <- ifelse(norw$r22meme_handshaking_community<97, norw$r22meme_handshaking_community, NA)

m.n <- lm(shake ~ as.factor(shakexp), data=norw, weights=norw$r22Weight2cap)

norw$shakexp.re  <- relevel(as.factor(norw$shakexp), ref=3)

m.n <- lm(shake ~ as.factor(shakexp.re), data=norw, weights=norw$r22Weight2cap)

w.norw <- svydesign(id      = ~responseid,
                    weights = ~r22Weight2cap,
                    data    = norw)

# Figure A2: Replication of the Handshaking Experiment using a probability sample from Norway. 

par(mfrow=c(1,1), mgp=c(2, .5, 0), family="Gill Sans", font.main=1, cex.main=1, cex.lab=1, cex.axis=1)
plot(c(1, 2),
     c(svymean(~shake, w.norw[norw$shakexp==3], na.rm=T)[1], svymean(~shake, w.norw[norw$shakexp==4], na.rm=T)[1]), pch=19, type="o", ylim=c(0, 1), xlim=c(0.7, 2.3), axes=F, ylab="Should Be Asked to Shake Hands [%]", xlab="", main="Replication: Norwegian Citizen Panel")

points(c(1, 2),
       c(svymean(~shake, w.norw[norw$shakexp==1], na.rm=T)[1], svymean(~shake, w.norw[norw$shakexp==2], na.rm=T)[1]), type="o")

abline(h=.5, col="lightgrey", lwd=1, lty=1)
axis(1, at=c(1, 2), labels=c("No Mention", "Hand on Heart"))
axis(2, at=c(0, .25, .5, .75, 1), label=c(0, 25, 50, 75, 100))

segments(1, svymean(~shake, w.norw[norw$shakexp==3], na.rm=T)[1]-1.96*SE(svymean(~shake, w.norw[norw$shakexp==3], na.rm=T)), 1, svymean(~shake, w.norw[norw$shakexp==3], na.rm=T)[1]+1.96*SE(svymean(~shake, w.norw[norw$shakexp==3], na.rm=T)))
segments(1, svymean(~shake, w.norw[norw$shakexp==1], na.rm=T)[1]-1.96*SE(svymean(~shake, w.norw[norw$shakexp==1], na.rm=T)), 1, svymean(~shake, w.norw[norw$shakexp==1], na.rm=T)[1]+1.96*SE(svymean(~shake, w.norw[norw$shakexp==1], na.rm=T)))
segments(2, svymean(~shake, w.norw[norw$shakexp==4], na.rm=T)[1]-1.96*SE(svymean(~shake, w.norw[norw$shakexp==4], na.rm=T)), 2, svymean(~shake, w.norw[norw$shakexp==4], na.rm=T)[1]+1.96*SE(svymean(~shake, w.norw[norw$shakexp==4], na.rm=T)))
segments(2, svymean(~shake, w.norw[norw$shakexp==2], na.rm=T)[1]-1.96*SE(svymean(~shake, w.norw[norw$shakexp==2], na.rm=T)), 2, svymean(~shake, w.norw[norw$shakexp==2], na.rm=T)[1]+1.96*SE(svymean(~shake, w.norw[norw$shakexp==2], na.rm=T)))


points(c(1, 2),
       c(svymean(~shake, w.norw[norw$shakexp==1], na.rm=T)[1], svymean(~shake, w.norw[norw$shakexp==2], na.rm=T)[1]), pch=19, col="white", cex=.7)
text(1, svymean(~shake, w.norw[norw$shakexp==3], na.rm=T)[1]-1.96*SE(svymean(~shake, w.norw[norw$shakexp==3], na.rm=T)), "Politicians", pos=1, cex=.8)
text(1, svymean(~shake, w.norw[norw$shakexp==1], na.rm=T)[1]+1.96*SE(svymean(~shake, w.norw[norw$shakexp==1], na.rm=T)), "Female Politicians", pos=3, cex=.8)

text(1, svymean(~shake, w.norw[norw$shakexp==3], na.rm=T)[1], round(svymean(~shake, w.norw[norw$shakexp==3], na.rm=T)[1]*100), pos=2, cex=.8)
text(1, svymean(~shake, w.norw[norw$shakexp==1], na.rm=T)[1], round(svymean(~shake, w.norw[norw$shakexp==1], na.rm=T)[1]*100), pos=2, cex=.8)
text(2, svymean(~shake, w.norw[norw$shakexp==4], na.rm=T)[1], paste(round(svymean(~shake, w.norw[norw$shakexp==4], na.rm=T)[1]*100), "**"), pos=4, cex=.8)
text(2, svymean(~shake, w.norw[norw$shakexp==2], na.rm=T)[1], paste(round(svymean(~shake, w.norw[norw$shakexp==2], na.rm=T)[1]*100), "**"), pos=4, cex=.8)
box()

##############################################################################################
### END
